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q ■ Abstract 

^ \ The Weber problem consists of finding a point in R n that minimizes the weighted sum of distances 
^ \ from m points in R n that are not collinear. An application that motivated this problem is the optimal 
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location of facilities in the 2-dimensional case. A classical method to solve the Weber problem, proposed 
by Weiszfeld in 1937, is based on a fixed point iteration. 

In this work a Weber problem constrained to a closed and convex set is considered. A Weiszfeld- 
like algorithm, well defined even when an iterate is a vertex, is presented. The iteration function Q 
\ that defines the proposed algorithm, is based mainly on an orthogonal projection over the feasible set, 
Q ■ combined with the iteration function of a modified Weiszfeld algorithm presented by Vardi and Zhang 
' in 2001. 

It can be seen that the proposed algorithm generates a sequence of feasible iterates that have descent 
properties. Under certain hypotheses, the limit of this sequence satisfies the KKT optimality conditions, 
is a fixed point of the iteration function that defines the algorithm, and is the solution of the constrained 
■ minimization problem. Numerical experiments confirmed the theoretical results. 
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1. Introduction 

Let a 1 , . . . , a m be m distinct points in the space R", called vertices, and positive numbers u>i, . . . , w m , 
called weights. The function / : R n — > R defined by 

m 

is called the Weber function, where ||-|| denotes the Euclidean norm. It is well-known that this function 
is not different iable at the vertices, and strictly convex if the vertices are not collinear (we will assume 
this hypothesis from now on). 

The Weber problem (also known as the Fermat- Weber problem) is to find a point in R n that min- 
imizes the weighted sum of Euclidean distances from the m given points, that is, we have to find the 
solution of the following unconstrained optimization problem: 

argmin f(x) 

x (2) 
subject to x G R n . 
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10 This problem has a unique solution x u in R n . 

11 The problem was also stated as a pure mathematical problem by Fermat [44, 27], Cavalieri [37], 

12 Steiner [14], Fasbender [20] and many others. Several solutions, based on geometrical arguments, were 

13 proposed by Torricelli and Simpson. In [30] historical details and geometric aspects were presented by 

14 Kupitz and Martini. In [41] Weber formulated the problem (2) from an economical point of view. The 

15 vertices represent customers or demands, the solution to the problem denotes the location of a new 

16 facility, and the weights are costs associated with the interactions between the new facility and the 

17 customers. 

is Among several schemes to solve the Weber location problem (see [12, 19, 28, 34]), one of the most 

19 popular methods was presented by Weiszfeld in [42, 43]. The Weiszfeld algorithm is an iterative method 

20 based on the first-order necessary conditions for a stationary point of the objective function. 

21 If we define T : R n ->■ R" by: 



if x 7^ a 1 , ... , a m , 

(3) 

a k , if x = a h , k = 1, . . . , m, 

22 the Weiszfeld algorithm is: 

x (l) = T , I G N, (4) 

23 where x^ G R n is a starting point. 

24 The Weiszfeld algorithm (4), despite of its simplicity, has a serious problem if some x^ lands ac- 

25 cidentally in a vertex a h , because the algorithm gets stuck at a k , even when a k is not the solution of 

26 (2). Many authors studied the set of initial points for which the sequence generated by the Weiszfeld 

27 algorithm yields in a vertex (see [29, 11, 6, 9, 7, 3]). Vardi and Zhang [40] derived a simple but nontrivial 

28 modification of the Weiszfeld algorithm in which they solved the problem of landing in a vertex. 

29 Generalizations and new techniques for the Fermat- Weber location problem have been developed in 

30 recent years. In [18] Eckhardt applied the Weiszfeld algorithm to generalized Weber problems in Banach 

31 spaces. An exact algorithm for a Weber problem with attraction and repulsion was presented by Chen et 

32 al. in [13]. Kaplan and Yang [24] proved a duality theorem which includes as special cases a great variety 

33 of choices of norms in the terms of the Fermat- Weber sum. In [10] Carrizosa et al. studied the so called 

34 Regional Weber Problem, which allows the demand not to be concentrated onto a finite set of points, 

35 but follows an arbitrary probability measure. In [17] Drezner and Wesolowsky studied the case where 

36 different l p norms are used for each demand point. In [23] the so called Complementary Problem (the 

37 Weber problem with one negative weight) was studied by Jalal and Krarup, and geometrical solutions 

38 were given. In [15] Drezner presented a Weiszfeld-like iterative procedure and convergence is proved if 

39 appropriate conditions hold. 

40 In some practical problems it is necessary to consider barriers (forbidden regions). Barriers were first 

41 introduced to location modeling by Katz and Cooper [25]. There exist several heuristic and iterative 

42 algorithms for single-facility location problems for distance computations in the presence of barriers (see 

43 [2, 8, 5, 4]). In [35] Pfeiffer and Klamroth presented a unified formulation for problems with barriers 

44 and network location problems. A complete reference to barriers in location problems can be found in 

45 [26]. Barriers can be applied to model real life problems where regions like lakes and mountains are 

46 forbidden. 

47 On the other hand, there are location problems whose solution needs to lie within a closed set. 

48 For example, see [39] for a discussion of the case when the solution is constrained to be within a 
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49 maximum distance of each demand point. Drezner and Wesolowsky [16] studied the problem of locating 

50 an obnoxious facility with rectangular distances (Zi norm), where the facility must lie within some 

51 prespecified region (linear constraints). A primal-dual algorithm to deal with the constrained Fermat- 

52 Weber problem using mixed norms was developed in [33] by Idrissi et al.. In [21] Hansen et al. presented 

53 an algorithm for solving the Weber problem when the set of feasible locations is the union of a finite 

54 number of convex polygons. In [36] Pilotta and Torres considered a Weber location problem with box 
constraints. 

56 Constrained Weber problems arise when we require that the solution is in an area (feasible region) 

57 determined by, for example, environmental and/or political reasons. It could be the case for a facil- 

58 ity producing dangerous materials that must be installed in a restricted (constrained) area. Another 

59 example could be the location of a plant in an industrial zone or of a hospital in a non-polluted area. 
In this paper a constrained location problem is considered. An algorithm is proposed to solve the 

61 following problem: 

argmin f(x) 

(5) 

subject to igO, 

62 where Q is a closed and convex set, generalizing the problem formulated in [36]. Problem (5) could be 

63 seen as a nonlinear programming problem and solved by standard solvers, but they may fail since the 

64 Weber function is not differentiable at the vertices. 

65 It can be proved that problem (5) has a unique solution x*, since the function / is strictly convex 
ee and Q is a closed and convex set. On the other hand, it is well-known that the convex hull of the given 

67 vertices a 1 ,. . . ,a m contains the solution x u of the unconstrained Weber problem (see for instance [29, 

68 pp. 100]). If Q contains the convex hull, both solutions x* and x u agree. In other cases, the solution 

69 x* is not necessarily a projection of x u over Q (see [36]). The algorithm is based basically on a slight 
variation of an orthogonal projection of the Weiszfeld algorithm presented in [40], that is well defined 
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71 even when an iterate coincides with a vertex. Properties of the sequence generated by the proposed 

72 algorithm related with the minimization problem 5 will be proved in the following sections. 

73 The paper is structured as follows: Section 2 describes the results in [40] in which a modified 

74 Weiszfeld algorithm is presented and some notation is introduced. In Section 3 the proposed algorithm 

75 is defined. Section 4 is dedicated to definitions and technical lemmas. In Section 5 the main results 

76 about convergence to optimality are presented. Numerical experiments are considered in Section 6. 

77 Finally, conclusions are given in Section 7. 

78 Some words about notation. As it was mentioned, we will call x u the solution of problem (2) and 

79 x* the solution of problem (5). The symbols || • || and (■, ■) will refer to the standard Euclidean norm 
so and standard inner product in R n , respectively. For a function / : R -> R we will denote by f'(a-) the 
si left-hand side derivative at a, and by f'(a+) the right-hand side derivative at a. 



82 2. The modified Weiszfeld algorithm 

83 This section reviews the main results presented in [40] in which the authors generalize the Weiszfeld 

84 algorithm for the case that an iterate lands on a vertex. From now on, this algorithm will be referred 

85 to as the modified Weiszfeld algorithm. 

se In order to make notation easier, we define the function A : R™ — > R by: 

f m 



A(x) 



> -n — 3 —rr, if x ^ a 1 , . . . ,a m , 
^1\\x-a?\Y T ' ' ' 



£ 



(6) 



ll x — a , k — 1, . . . , m. 



2\\a k - a'J 
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Notice that A(x) > for all x G R™. In [40, pp. 563], the number A(a k ) was called A^. 

A generalization for the iteration function T , defined in (3), is given by T : R n — > R™ defined as 
89 follows: 

Wjd 1 



f m 



T(x) 



E 



\x — a ] \ 



2A(x) 



if x 7^ a 1 



E 



.„3 



(7) 



a fe — (P\ 



2A (a k 



if x = a h , k = 1, 



m. 



90 Notice that T coincides with T in R n — {a 1 , . . . , a m }. 
si Let R : R n R n and r : E n -> E be: 



E 

j'=i 

m 

E 



(a- 7 — x) 
|x — a- 7 II 



(a J — a fc ) 



if x 7^ a 1 , . . . , a r 



if x = a , k = 1, . . . , m, 



r(x) 



\\R(x)\\, VxGR™. 

92 The function R generalizes the negative gradient of the Weber function since, for all x 7^ a 1 

V/(x) = -R(x). 



m 

)•••)« j 



(9) 



93 The following lemma is very easy to prove (see [40, equation (14)]), and it relates the functionals T 

94 and R. 

95 Lemma 1. For all x G R n we have R(x) = 2A(x) T(x) — x 

If we define 7 : R n — > R by: 

0, if x 7^ a 1 , ... , a m , 

j(x) = \ 0, if x = a fc and r (a fc ) = for some k = 1 

■y^/ 7 " , if x = a k and r (a fc ) 7^ for some k = 1 



96 we can see that 7(2) > for all x G R n . 

The modified Weiszfeld algorithm presented in [40] is defined by: 

x (D=T(x^), I G N, 

97 where x^ G R n and T : R n — > R n is given by: 

T(x) = (1 - /3(x)) f (x) + /3(x)x, 

98 where (3 : R™ — >■ R is defined by (3(x) = min {1, 7(x)}. 



(10) 



99 Remark 2. (aj If x ^ a 1 ,. . . ,a m , then /3(x) = because 7(x) = 0. So, we can deduce that T(x) 

100 T(x). Notice that this fact implies that the functional T is continuous in R n — {a 1 , . . . , a m }. 
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(b) It can be seen that if a k 7^ x u , then < (3(a k ) < 1 (see [40, pp. 563]). 

(c) From equation (10) we obtain that T(x) — x = (1 — f3(x)) (t(x) — x\ for x G R". 



103 The main result in [40, pp. 562] is: 

104 Theorem 3. The following propositions are equivalent: 

105 (a) x = x u . 

106 (b) T(x) = x. 

107 (c) r(x) < rj(x). 

where 

0, if x 7^ a 1 , . . . ,a r 



r](x) 



Wk, if x = a , k = 1, . . . , m. 



108 



3. The proposed algorithm 



109 This section is dedicated to describe the proposed algorithm, introducing some definitions and re- 
no marks. 

in First of all, we can notice that problem (5) has a unique solution, due to the fact that / is a 

112 non-negative, strictly convex, and continuous function, lim|| x ||_> 00 f(x) = 00 and Q is closed and convex. 

In order to define the proposed algorithm at the vertices, we will need to determine which points of 
the segment that joins a k and T(a h ) are in the feasible set Q. If k = 1, ... , m, let the set Sk be defined 
by: 

S k = {X G [0, 1] : (1 - X)T(a k ) + Xa k G fi} . 

Notice that Sk could be equal to the empty set in case that a h and T(a k ) do not belong to Q. On the 
other hand, if a k G f2, then 1 G Sk, which means that Sk 7^ 0. Thus, we can define: 

X(a h ) = inf S k , a h G fi. 

113 In case a vertex a k is not in Q, there is no need to define the number X(a k ). 

114 In the following lemma, a set of basic properties of X(a k ) are listed: 

115 Lemma 4. If k = 1, . . . ,m and a k G Q then: 

lie (a) X(a k ) G [0,1]. 

117 (b) IfT(a k ) G tt then X(a k ) = 0. 

us (c) IfT(a k ) £ Q then X(a k ) G (0, 1]. 

119 Proof. The proof of (a) follows from the definition of S k . If T(a k ) G Vt, then G S k , so X(a k ) = 0, 

120 and this proves (b). Finally, for item (c), let us consider that T(a k ) ^ Q. Since Q is a closed set, there 

121 is an entire ball centered at T(a k ) that does not intersect Q, which implies that there exists e such that 

122 (1 - X)T(a k ) + Xa k £ for all A G [0, e]. Thus, X(a k ) G (0, 1] and this concludes the proof. □ 

123 Let us call Pq : R™ — > Q the orthogonal projection over Q. Since Q is a nonempty, closed and convex 

124 set, the operator Pq is a continuous function [1, pp. 99]. 

125 We define the iteration function Q : Q — > Q by: 

C>(r\ - S p n° T (x), iix^ a l ,...,a m , , , 

^ { ) \ (1- X{a k ))T(a k ) + X(a k )a k , if x = a k G ft, k = 1, . . . , m. 1 ' 
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126 There will be no need to define Q outside f2 since the proposed algorithm generates a sequence of feasible 

127 points. The iteration function Q at x G f2 coincides with the orthogonal projection of T(x) over the 

128 feasible set when x is different from the vertices. Only when x is a vertex a k belonging to Q, Q(x) is 

129 defined as the farthest possible feasible point of the segment that joins x with T(x). 

130 The following remark states some basic properties of the iteration function of the proposed algorithm. 



131 



132 Remark 5. 

133 (a) Ifa k G Q and T(a k ) G Q, then Q(a k ) = T(a k ) =P n o T(a k ). 

134 (b) If a k G Q, it can be seen that: 

Q(a k ) -a k = (1 - \(a k )) (T(a k ) - a k ) , 
Q(a k ) - T(a k ) = -X(a k ) (T(a k ) - a k ) . 

135 (c) The functional Q is continuous in IR™ — {a 1 , . . . , a m }. 

136 Proof. The proofs of (a) and (b) are straightforward. For (c), since Pq, is continuous in R n (see [1, 

137 pp. 99]) and T is continuous in R" — {a 1 , . . . ,a m } (see Remark 2), we have that Q is continuous in 

138 R n - {a\...,a m }. □ 

139 The proposed algorithm is described below. 

wo Algorithm 6. Let f2 C R, n be a closed and convex set. Assume that x (0) G tt is an initial approximation 

141 such that f(x^ ') < /(a J ) for all j G {1, . . . ,m} and a? G f2. Given e > a tolerance and x^ -1 - 1 G f2 ; do 

142 the following steps to compute x^ : 

143 

144 Step 1: Compute: 

Step 2: Stop the execution if 



« = Q . (12) 



x 



\ X V)- X «-V\\ < e , 



145 and declare x^ 1 ' as solution to the problem (5). Otherwise return to Step 1. 

146 From the definition of Q it follows that Algorithm 6 generates a sequence of feasible iterates. Also 

147 notice that if there are vertices in the feasible set, x^ can be one of them, for example, a vertex a s such 

148 that f(a s ) < for all a? G f2. On the other hand, if there are no vertices in the feasible set, x^ can 

149 be chosen as the projection over Q of the null vector. 

150 4. Some definitions and technical results 

151 The purpose of this section is to define some entities and prove technical lemmas that will be 

152 important in the proof of the main results. 

First of all, we will define some useful operators for making notation easier. If A C {1, . . . , n}, then 
we define || • \\ A : R n ->■ E and (•, -) A : R n x R" -> R by: 

153 Notice that, when A C {1, . . . , n}, || • \\ A is not necessarily a norm and (-, -) A is not necessarily an inner 

154 product. 
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155 According to this definition, if A and B are sets such that A D B = and AU B — {1, . . . , n}, it can 
we be seen that: 



INI 

c(x,v)a 



\ x \\a + Hulls' 



(cx,y) A = (x,cy) A . 



(13) 
(14) 
(15) 



157 For x G O, let us define the following sets of indices: 

Af(x) = {fceN: l<k<n, (T(x)) k ^(Q(x)) k }, 
S(x) = {fceN: l<A;<n, (T(z)) fc = (Q(z)) fc } , 

15 8 Notice that for all x G R n we have that A/"(x) n £(ar) = and Af(x) U£(x) = {1, 

159 Let a : f2 — > R" be the following function: 



n}. 



If x ^ a 1 , 



x — a-? 



[g(x) - o*] . 



(16) 



If x = a k G Q for some k = 1, . . . , m: 



, . . . , 

m 



j'=i 



a" — a J 



/3(a fc ))a J - f3{a k )a k ] . 



(17) 



162 It can be seen that the function a is related to the iteration function Q of the proposed algorithm, 

163 and the iteration function T of the modified algorithm. 



164 Lemma 7. If x G ft, then a(x) = 2A(x) [Q{x] - T{x)\ . 

165 Proof. If x ^ a 1 , . . . , a m , then: 



a(x) 



E 



x — a J 



E 



\x — a ] \ 



g(x) 



i — a- 7 ^— f kr — a J 

Em 



x— a' 



2A(x) [Q(x) - T(x)] 



Em 10 



x—a,3 



2A(x) Q(x) - T(x) 



lee where in the last equalities we have used the definition of T as in (7), and the fact that T(x) = T(x) 
167 due to Remark 2. 

lea If x = a k for some k — 1, . . . , m, we follow a similar procedure than in the previous case. □ 



169 
170 



Now, we will define auxiliary functions that take into account the projection Pq in order to prove a 
descent property of / (see next sections). If x G Q, we define: 



i7i (a) E x : R n -> R", where: 



(Q(x)) fc , iikeAf(x), 
y k , iikeS(x). 
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(b) If £(x) = {ix, . . . , z r } ^ define P x : R n ->■ R r where: 

( p Ay))k = Vik> k = l,...,r. 
A useful property of E x , that follows from the definition, is pointed out in the following remark. 
Remark 8. If x G then E x o Q(x) = Q(x). 

The iteration function Q inherits an important property from the orthogonal projection Pq. 
Lemma 9. If x G Q we have that (Q(x) — x, Q(x) — T(x)) < 0. 

Proof. If x ^ a , ... , a m , then Q(x) = Pq o T(x). By a property of the orthogonal projection [1, pp. 
93] we have that (Q(x) - x, Q(x) - T{x)) < 0. 

If x = a k for some k — 1, . . . , m, Remark 5 and Lemma 4 imply: 

(Q(a k ) - a\ Q(a k ) - T(a k )) = -X(a k ) (l - X(a k )) \\T(a k ) - a k \\ 2 < 0, 
and this concludes the proof. □ 

The next technical lemma will help us to save computations in other lemmas. 
Lemma 10. If x G f2, A is a subset of {1, ... ,n} and j G {!,..., rn}, then: 



\\Q(x) - a 3 \\ A 
Proof. If x G Q, we have: 
\\Q(x)-ai\\ 2 A 



\x — a J 



- \\Q{x) - x\\\ + 2 (Q(x) - x, Q(x) - a j ) 



A 



(Q(x) - a J ,Q(x) - a ] ) A 

(Q(x) — x + x — a? , Q(x) — x + x — ai) A 



\\Q(x) — x\\ A + ||x — a? ||^ + 2(Q(x) — x,x — a J ) 



A 



||Q(x) -xf A + \\x-a j \\ 2 A + 2(Q{x) -x,x-Q{x)) 



A 



+ 2{Q(x) -x,Q(x) -a j ) 



A 



\x — a 



in 2 
U 



\\Q{x) -x\\ 2 A + 2{Q(x) - x, Q(x) - a j ) 



A- 



iB2 If x G Q, let us define g x : R n —> R by: 



9x(y) = < 



E 

3+k 



.ill 2 



if x 7^ a 1 , . . . , a m 



2\\a k - aP\ 



\y — a J || 2 + Wk\\y — a k \\, if x = a k , 



k = 1, . . . , m. 



□ 



(19) 



183 The values that g x assumes at x and Q(x) will play an important role in the proof of a property of the 

184 objective function /. 



185 Lemma 11. Let x G Q be. 



186 (a) If x ^ a 1 , . . . ,a m then: 

g x (x) = ^f{x) + 2A{x){Q{x)-x,Q(x)-T(x)) 
- A{x)\\Q{x)-x\\ 2 N{x) . 

187 (b) If x = a k for some k — 1, . . . , m, then g a k(a k ) = \f(a k ). 

Proof. Let us suppose that i/a 1 ,..., a m . By property (13) and (18), we have for j — 1, . . . , m: 



\E x (x) - a j \\' = \\x - a ] \\ 2 £(x) + \\Q{x) - a j \\^ {x) 



Using Lemma 10, we can see that: 

m 

.7=1 11 1 



I 3 1 1 ^ _i_ 1 1 3 I I ^ 

r ~~ a \\s(x) + \\ x - a \\m( x ) 



- \\Q{x) - x\\ 2 N{x) + 2 (Q(x) - x, Q(x) - a j ) M{x) 

Due to (13), the definition of the Weber function /, the definition of A as in (6), the property (15) and 
the definition of a as in (16), we obtain: 

9x{x) = -f(x) -A(x)\\Q(x) -x\\% {x) + (Q(x) - x,a{x)) M{x) . 

189 By Lemma 7, the fact that (Q(x)) i = (T(x))^ for all i G £{x) and (14), we get: 

9x{x) = ^f(x) - A{x)\\Q(x) -xf N{x) 
+ 2A(x)(Q(x)-x,Q(x)-T(x)), 

190 which concludes the proof of (a). 

Now, let us assume that x = a k for some k — 1, . . . , m. Then: 

9a^a k ) = ^_^-_||a fe -^|| 2 = -J2w j \\a k -a>\\ = -f(a k ). 

3=1 J=l 

191 This concludes the proof of (b). □ 

192 The number g x (Q(x)) can be computed in the next lemma. 

193 Lemma 12. Let x £ Q be. 

194 (a) If x ^ a 1 , . . . ,a m then: 

g x (Q(x)) = ^f(x)+2A(x)(Q(x)-x,Q(x)-T(x)) 
- A(x) \\Q(x) - xf . 

195 (b) If x = a k for some k — 1, . . . , m, then 



9a*(Q(a k )) = l -f{a k )-A{a k )\\Q{a k )-a k \\ 2 

+ 2A(a k )(Q(a k )-a k ,Q(a k )-T(a k )) 

- 2f3(a k )A(a k ){Q{a k )-a k ,f(a k )-a k 

+ w k \\Q{a k ) - a k \\ . 
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Proof. First, let us consider x ^ a 1 , . . . , a m . Due to Remark 8 we have: 



g x (Q(x)) = p^ ^^ \\Q(x)-ai\\ 2 . 



196 By Lemma 10 we obtain: 



9*(Q(*)) = E on W ' J 3 -n ^-* j \\ 2 -\\Q{x)-A? 

' 2 x — a 3 \\ 
j=i ii 

+ 2(Q(x) -x,Q(x) -a J >] . 

Due to the definition of the Weber function /, the definition of A as in (6) and the definition of a as in 
(16), we deduce that: 

9x(Q(x)) = \f{x) - A(x)\\Q{x) -x\\ 2 + (Q(x) - x, a{x)). 

By Lemma 7 we get: 

9x(Q(x)) = \f{x) - A(x)\\Q(x) - x\\ 2 + 2A(x)(Q(x) - x, Q(x) - T(x)), 

197 concluding the proof of (a). 

Now, consider x = a k for some k — 1, . . . , m. Due to (19) we have: 

m 

9 AQ(A) = E on P J IQfo*) - aJ W 2 + MQ^ k ) -AY 

z — ' 2 kr — a? 

j=i " " 

a* 

198 By Lemma 10, the definition of the Weber function / and the definition of A as in (6) we obtain: 

9 ak {Q{a k )) = l -f(a k )-A{a k )\\Q{a k )-a k \\ 2 

I m 

+ w k \\Q(a k ) - a k \\. 

Manipulating algebraically, 

Q(a k ) - a J = Q(a k ) - (1 - (3(a k ))a J - (3(a k )a k + (3(a k )(a k - a ] ). 

199 Due to the definition of a (see (17)) and the definition of R (see (8)) we get: 

9a*{Q(a k )) = l -f(a k )-A(a k )\\Q(a k )-a k \\ 2 

+ (Q{a k ) - a k , a{a k )) - f3(a k ) (Q(a k ) - a k , R{a k ) 
+ w k \\Q(a k ) -a k \\. 

200 By Lemma 1 and Lemma 7 we have: 

9a^{Q(a k )) = l -fXa k )-A(a k )\\Q(a k )-a k \\ 2 

+ 2A(a k ) (Q(a k ) - a k , Q(a k ) - T(a k )) 

- 2A(a k )(3(a k ) (Q{a k ) - a k ,f(a k ) - a k ^ + w k \\Q(a k ) - a k \\. 

201 which concludes the proof. □ 
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202 The next lemma deals with the last two terms of g a k {Q{ ak ))- 
Lemma 13. If a k G Vt for some k = 1, . . . , m, the number 

z = Wk \\Q(a k ) - a k \\ - 2A(a k )(3(a k ) (Q{a k ) - a k ,f(a k ) - c 

203 is equal to zero. 

204 Proof. If Q(a k ) = a k the result is true. So, from now on, let us consider that Q(a k ) ^ a k . First, let us 

205 check that a k ^ x u . In case that a k = x u , then T(a k ) = a k by Theorem 3. Since a k G f2, then T(a k ) G f2. 

206 By Remark 5 we have that Q{a k ) = T(a k ) = a k which is a contradiction. 

By Remark 2, we have that f3(a k ) G (0, 1) (since a k ^ x u ) and: 



w k \\Q(a k 



o*|| _ 2A ( ak )P( ak ) ( Q{ak) _ akj T{ak) _ 



1 - (3{a k ) 



207 Extracting common factors, using Remarks 2 and 5, the fact that T(a k ) ^ a k (if T(a k ) = a k then 

208 a k = x u by Theorem 3), and the fact that f(a k ) ^ a k (if f(a k ) = a k then T(a k ) = a k by definition (10)) 

209 we get that: 



2A(a k ) \\Q(a k 



T(a k 



2A(a k 



T(a k ) 



(3(a k ) 



(1 - X(a k ))(T(a k ) -a k ) T(a k ) - a k 



(1 - X(a k )){T{a k ) - a k )\\ ' \\T{a k ) - a k \ 
Simplifying and using the definition of (3(a k ) we have that: 



z = 2A(a k ) \\Q(a k ) 
210 which concludes the proof. 



T{a k ) - a k [(3(a k ) - f3{a k )\ = 0, 



□ 



211 The purpose of the next two lemmas is to determine a strict inequality between the functions g x and 

212 / at suitable points. First of all, we have to prove the following result. 

213 Lemma 14. Let x G Q be such that x ^ Q(x). 

214 (a) If x ^ a 1 ,. . . ,a m , then g x (Q(x)) < 9x(%)- Besides that, if £{x) ^ and P x o Q(x) ^ P x (x), then 

215 9x(Q(x)) < 9x(x). 

216 (b) If x = a k for some k = 1, ... ,m, then g a k(Q(a k )) < g a k(a h ) . 



217 Proof. If x ^ a l ,...,a m , then g x (Q(x)) — g x {x) 



-A(x)\\Q(x) 



x 



£(x) 



< 0, by Lemma 11 and 



218 Lemma 12. Besides that, if S(x) ^ and P x o Q(x) ^ P x {x) we deduce that ||Q(x) — x||£( x ) ^ 0. Thus, 

219 g x (Q(x)) < g x (x). 

220 If x = a k for some k — 1, . . . , m, by Lemmas 11,12 and 13 we have: 

g a *(Q(a k )) - g ak (a k ) = -A(a k )\\Q(a k )-a k \\ 2 

+ 2A(a k )(Q(a k )-a k ,Q(a k )-T(a k )). 

Due to Lemma 9 and the fact that A > we obtain: 

g a *(Q(a k ))-g ak (a k ) < -A(a k )\\Q(a k ) - a k \\ 2 < 0, 

221 and the proof is finished. □ 
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222 Lemma 15. Let x G Q be such that x 7^ Q(x). Then g x (Q(x)) < ^f(x). 

223 Proof. Let us consider the case when i^a 1 ,..., a m . By Lemmas 9, 11 and 14 we have that: 

9x{Q{x)) < g x {x) = ^f(x) + 2A{x){Q{x)-x,Q{x)-T(x)) 
- A(x) \\Q(x) -xf M{x) 
< ±f(x)-A{x)\\Q(x)-x\\*f ilB y 
If £{x) = 0, then || • \\u{x) — II ■ II- Therefore: 

9*(Q(x)) < \f{x) - A(x)\\Q(x) -x|| 2 < l/(x). 

If £(2) 7^ and P x o = P x (x), then there exists an index % G jV(:r) such that Xj 7^ (Q(a;))j since 
x 7^ Thus, ||Q(sc) — ^||jV(x) 7^ 0, which implies: 

g x (Q(x)) < i/(x) - ||Q(x) - zll^j < -/(*). 

If 7^ and P x o 7^ P x (x), due to Lemmas 9, 11 and 14, we have that: 

g x {Q{x)) < g x (x) < l -f{x) - A(x) \\Q(x) - x\\ 2 M[x) < l -f{x). 

224 Now, when x = a k for some k = l,...,m, g a k(Q(a k )) < g a k(a k ) = |/(a fc ) due to Lemma 11 and 

225 Lemma 14. □ 

226 The next lemma states an equality that relates the Weber function and g x at appropriate points 

227 when x 7^ a 1 , . . . , a m . Besides that, this result will be crucial in the next section. 

Lemma 16. Let x 7^ a , . . . , a m be such that x G Q and x 7^ Q(x). Then: 

g x ° Q(x) = + (f(Q(x)) - f(x)) + 5, 5 > 0. 

Proof. Due to the definition of g x as in (19) and Remark 8 we get that: 

m 



g x o Q{x) = V — — 3 -—\\Q(x) - a 
*— ' 2 1 — a 3 \\ 



Adding and subtracting \\x — a J \\ we have: 



9x°Q(x) = ^^rr^^[\\x ~a J \\ + (\\Q{x) -a 3 \\-\\x -a 3 ^)]' 



^ m m 
3=1 3=1 



a J \\ — \\x — a 3 



j=l 



Notice that the first term of the last equality is the Weber function (divided by two), and the last term 
is a non-negative number, so we will define it as 5. So, using the definition of the Weber function in the 
middle term we obtain: 

g* ° Q(x) = \f{x) + (f(Q(x)) - f(x)) + 5. 

□ 
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230 5. Convergence to optimality results 

231 This section states the main results about convergence of the sequence {a:®} generated by Algorithm 

232 6. The next theorem establishes that if a point x 6 VL is not a fixed point of the iteration function, then 

233 the function / strictly decreases at the next iterate. 

234 Theorem 17. Let x E f2 be such that x ^ Q(x). Then f(Q(x)) < f(x). 
Proof. Let us consider that x ^ a 1 , . . . , a m . By Lemma 15, we have that: 

g x oQ( x ) < ^f(x). 

By Lemma 16 we get that: 

lf(x) + f(Q(x))-f(x) + 5<^f(x). 
Simplifying the last expression we obtain: 

f(Q(x))-f(x) + 5<0. 

Finally, 

f(Q(x)) - f(x) < f(Q(x)) - f(x) + 5<0. 

235 Therefore, f(Q(x)) < f(x). 

Now, consider that x = a k for some k = 1, . . . ,m. Following a reasoning similar than in [40, pp. 
564], using Lemma 14 we have that: 

g ak oQ{a k ) -g ak (a k ) < 0. 

236 By definition of g a k we know that: 

g ak oQ( a k ) -g a k(a k ) = w k \\Q(a k ) - a k \\ 

m 

+ E^^(n^)-^ir-iK-^ir)- 

j=i " " 

237 Using the fact that (a 2 — b 2 ) / (2b) > a — b for a = \\Q(a k ) — a J || 2 > and b = \\a k — a j \\ 2 > we obtain 

238 that: 

g ak o Q( a k ) - g a k(a k ) > w k \\Q(a k ) - a k \\ 

m m 

- \\a k - a j \\ + ||(5(a fc ) - a j \\ . 

3=1 3=1 

239 Rearranging terms we deduce that: 

> g ak oQ(a k ) - g ak (a k ) 

m m 

= \\Q{a k ) - a j \\ - \\a k - a j \\ = f(Q(a k )) - f{a k ), 

3=1 3=1 

240 and the proof is complete. □ 
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241 Corollary 18. Let jx^} be the sequence generated by Algorithm 6. Then the sequence {/(^ )} is 

242 not increasing. Even more, each time x"' 7^ Q (x"n the sequence strictly decreases at the next iterate. 

243 If the sequence {x^} generated by Algorithm 6 were not bounded, then we could choose a subse- 

244 quence {y^} such that — > 00. But this implies that f(y®) — >■ 00, which is a contradiction since the 

245 sequence {/ (2?™)} is not increasing. 

246 Remark 19. The sequence |x®} generated by Algorithm 6 is bounded. So, there exists a subsequence 

247 convergent to a point x* G Q. Hence, x* is a feasible point. 

Due to the nondifferentiability of / at the vertices a 1 ,. . . , a m , we can not use the KKT optimality 
conditions at a k . Therefore, if a k and z are in Q, let us define G z k : [0, 1] — Y R by: 

G z ak {t) = f{a k + t{z-a k )). 

248 If a G f2, z G f2, t G [0, 1] and Q convex, we have that a k + t(z — a k ) G f2. Notice that the right-hand 

249 side derivative G^ k (0+) (or the directional derivative of / in the direction of z) exists (see [22, pp. 33]). 

250 Besides that, 



G z J(0+) =w k \\z-a k 



R(a k ), z- a k 



(20) 



251 The next lemma shows that if we are in a vertex a k , the directional derivative of / at a k in the 

252 direction of Q(a k ) is a descent direction. 



253 Lemma 20. Let a k G fl be such that T(a k ) Q. Then: 

G* ak '(0+) > G Q J ak) '(0+), Vz G [a k ,T(a k )} 

254 where: 



G 



Q(a k ) 1 



(0+) = -2 [1 - (3(a k )] A{a k ) T(a k ) 



\Q(a k ) 



(21) 



(22) 



255 Proof. If T(a k ) = a k then T(a k ) G Q, which is a contradiction. Besides that, if T(a k ) = a k , we would 

256 have that T(a k ) = a k because of (10), and again it would be a contradiction. So, we will consider 

257 T(a k ) ^ a k and T(a k ) 7^ a k for the rest of the proof. Since T(a k ) 7^ a k , then (3(a k ) G (0, 1) (see Remark 

258 2 and Theorem 3). 

Let us prove equation (22) first. Now, by (20) we can see that: 

G Q J ak) '(0+) = w k \\Q(a k ) - a k \\ - (R(a k ),Q(a k ) - a k 



Notice that if Q{a k ) = a k , equation (22) holds. So, let us consider from now on that Q{a k ) 7^ a k . By 
using Lemma 1 we replace R(a k ) and get: 



Q{a k ) I 



(0+) = w k \\Q{a k ) - a k \\ - 2A{a k ) (T{a k ) - a k , Q{a k ) 



G 



259 Extracting common factors and using the definition of (3 when it belongs to (0, 1) we obtain: 



G 



Q(a k ) I 



(0+) = 2A{a k ) \\Q(a k ) - a k \\ T(a k ) - a 



(3(a k ) 



T{a k 



Q{a k 



T(a k ) — a k 



' \\Q(a k ) ~a k \ 
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By Remarks 2 and 5 the vectors Q{a k ) — a k and T(a k ) — a k are parallel, so: 

G Q(a") / (Q + ) = 2A(flfc) ||g (afc) _ a fc|| _ a k J^^fcj _ X j _ 

260 which is equivalent to (22). 

Now, let us prove (21). If z = a k then G z ak {t) = f(a k ) for all t G [0,1], thus G*^(0) = 0, and 
therefore the inequality (21) holds. So, let us assume that z ^ a k for the rest of the proof. Using (20) 
and due to Lemma 1 to replace R(a k ): 



G z J(0+) =w k \\z- a k \\ - 2A(a k ) (T(a k ) -a k ,z- a 1 



26i Extracting common factors: 



G z J{0+) = 2A(a k )\\z-a k \\ T{a k ) - a k 



T(a k ) — a k z — a k 



2A(a k ) T(a k ) - a k 



T(a k ) — a k 



\z — a* 



262 Using the expression for (3(a k ) G (0, 1) we obtain: 



G* k '(0+) = 2A(a k )\\z-a k \\ T{a k ) - a 



T(a k ) -a k z-a h 



(5{a k ) 



T(a k ) - a h 



\z — a h 



If z belongs to the segment that joins a k and T(a k ) we have that z — a k and T(a k ) — a k are parallel 
vectors, then: 

G z J(0+) > -2 [1 -(3(a k )] A(a k ) \\z - a k \\ f(a k ) - a k . 
We can write z — (1 — X)T(a k ) + Xa k where A G [0, 1]. Therefore: 

G z J{0+) > -2 [1 - (3(a k )] A(a k ) (1 - A) \\T(a k ) - a k \\ f(a k ) - 

for all A G [0,1]. The minimum value of the right-hand side of the last expression happens when 
A = A(a fc ), so: 

G z J(0+) > -2 [1 - f3(a k )] A(a k ) (l - X(a k )) \\T(a k ) - a k \\ f(a k ) 
Using Remark 5 we conclude that: 

G z J(0+) > -2 [1 - (3(a k )] A(a k ) \\Q(a k ) - a k \\ f(a k ) - a 



□ 
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264 Now we will prove an equivalence that characterizes the solution of (5) in terms of the iteration 

265 function Q. Moreover, if x* is a regular point that is not a vertex, then x* is a KKT point. 

266 From now on, let us consider that 

n = {yeR n :g(y)<0,h(y) = 0}, (23) 

267 where g : R n — > R s is a convex function and h : R n — > W is an affine function. 

268 Theorem 21. Let Q be defined as in (23) and x G VL. Consider the following propositions: 

269 (a) x is a KKT point. 

270 (b) x is the minimizer of the problem (5). 

271 (c) Q(x) = x. 

272 If x 7^ a 1 , . . . ,a m , g and h are continuously differentiable, and x is a regular point, then (a), (b) and 

273 (c) are equivalent. 

274 If x = a k for some k = 1, . . . , m, then (b) implies (c). 

275 Proof. Let x ^ a 1 ,...,a m be. Since / is strictly convex and Q is convex, the KKT optimality 

276 conditions are necessary and sufficient. Therefore, it holds that (a) is equivalent to (b). 

277 Now we will prove that (b) implies (c). Let us suppose that x is the minimizer of the problem (5). 

278 If x were not a fixed point of the iteration function Q, we would have that x ^ Q{x), which means that 

279 f(Q( x )) < f( x ) by Theorem 17. This contradicts the hypothesis. 

To demonstrate that (c) implies (a), we will assume that x is a fixed point of Q, that is, x = Q(x). 
Since Q(x) = Pn ° T(x), x is the solution of: 

argmin F(z) = -\\z — T(x)|| 2 
z 2 

subject to g(z) < 0, 
h(z) = 0. 

280 Since F and g are convex, h is affine, and x is a regular point, the KKT optimality conditions hold at 

281 x. That is, there exist multipliers {/ij}* =1 and {Aj}J =1 such that (see [31, 38]): 

s p 

x — T(x) + fijVgj(x) + \jVhj(x) = 0, 

3=1 3=1 

Vj9j(x) = 0, j = l,...,s, 
Hi > 0, j = l,...,s, 

g(x) < o, 

h{x) = 0. 

282 Multiplying these equations by 2A(x), using equation (9), Lemma 1 and Remark 2, we obtain: 

s p 

V/(x) + ^(2A(x)/i J )V^(x) + ^(2A(x)A J )V/ ij (x) = 0, 

3=1 3=1 

(2A(x)iJ lj )g j (x) = 0, j = l,...,s, 
{2A{x) N ) > 0, j = l,...,s, 
g(x) < 0, 
h(x) = 0. 

283 where {2A(x)[ij} s =1 and {2A(x)Xj} p =1 are multipliers. Therefore, x is a KKT point of the problem (5) 

284 (see [31, 38]). 

285 Now, let us suppose that x = a h for some k — 1, . . . , m. As before, if x is a minimizer of the problem 

286 ( 5), then Q(a k ) = a k , otherwise f(Q(a k )) < f(a k ), which would be a contradiction. □ 
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287 6. Numerical experiments. 



2B8 The purpose of this section is to discuss the efficiency and robustness of the proposed algorithm 

289 versus a solver for nonlinear programming problems. 

290 A prototype code of Algorithm 6 was programmed in MATLAB (version R2011a) and executed in a 

291 PC running Linux OS, Intel(R) Core(TM) i7 CPU Q720, 1.60GHz. 

We have considered a closed and convex set Q C R 2 defined by the set Q = {y G R n : g(y) < 0}, 
where g is given by: 



9{x) 



17 1 

-4 x H x 2 H x 2 (x — 3) + y 

8 72 216 V 7 y 



-x 



59 
10 



11 

T 



3 

2 X 



35 

T 

13 



x-y- 



-4 + -(x - 1) + — (x - l) 2 + — {x - l) 2 (x - 3) - y 
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1 

~3 X 
2 

— x 
3 

-4x 



32 
11 

V ~T 

13 

V 

y 3 

y -19 



292 The feasible set is defined by linear and nonlinear constraints, as it can be seen in Figure 1. 



293 



294 We have built 1000 different experiments where for each one: 

295 • The number of vertices was m = 50. 

296 • The vertices were normally distributed random vectors, with mean equal to and standard devi- 

297 at ion equal to 10. 

298 • The weights were uniformly distributed random positive numbers between and 10. 

299 • Tolerance was set to e — 0.00001. 

300 On one hand, each experiment was solved using Algorithm 6 and, on the other hand, it was considered 

301 as a nonlinear programming problem and solved using function fmincon (see [32] and references therein) . 

302 Since the Weber function (1) is not different iable at the vertices, nonlinear programming solvers may 

303 fail. 

304 Let x m (i) be the solution of (5) obtained by fmincon in experiment i, and f m (i) = f(x m (i)). 

305 Analogously, let x p {i) be the solution of (5) obtained by Algorithm 6 in experiment i, and f p (i) = 

306 f(x p (i)). Figure 2 shows the difference between the arrays f m and f p . Both methods finished succesfully 

307 in all cases, however, Algorithm 6 found equal or better results for all experiments. For example, 

308 the difference f m — f p was greater than 0.01 in 35 experiments (the maximum difference ocurred in 

309 experiment 506). 
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Figure 1: Feasible set f2 
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Figure 2: Difference between minimum values found by Algorithm 6 and fmincon. 
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Feasibility of the solution 

0.5 r 
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Figure 3: Feasibility of the solution x p (i) obtained by Algorithm 6. 

310 

311 Feasibility of the solutions x p (i) can be checked computing m&x(g(x m (i))). Results can be seen in 

312 Figure 3 

313 

314 7. Conclusions 

315 This paper proposes a Weiszfeld-like algorithm for solving the Weber problem constrained to a closed 

316 and convex set, and it is well defined even when an iterate is a vertex. The algorithm consists of two 

317 stages: first, iterate using the fixed point modified Weiszfeld iteration (10), and second, either project 

318 onto the set Q when the iterate is different from the vertices, or, if the iterate is a vertex a k , take the 

319 point belonging to the line that joins T(a k ) with a k as defined in (11). 

320 It is proved that the constrained problem (5) has a unique solution. Besides that, the definition 

321 of the iteration function Q allows us to demonstrate that the proposed algorithm produces a sequence 

322 {a;®} of feasible iterates. Moreover, the sequence {/ (x^)} is not increasing, and when ^ Q (x^), 

323 the sequence decreases at the next iterate. It can be seen that if a point x* is the solution of the problem 

324 (5) then x* is a fixed point of the iteration function Q. Even more, if x* is different from the vertices, 

325 the fact of being x* a fixed point of Q is equivalent to the fact that x* satisfies the KKT optimality 

326 conditions, and equivalent to the fact that x* is the solution of the problem (5). These properties allows 

327 us to connect the proposed algorithm with the minimization problem. 

328 Numerical experiments showed that the proposed algorithm found equal or better solutions than a 

329 well-known standard solver, in a practical example with 1000 random choices of vertices and weights. 

330 That is due to the fact that the proposed algorithm does not use of the existence of derivatives at the 

331 vertices, because the Weber function is not different iable at the vertices. 
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